{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy import stats\n",
    "import math\n",
    "from data import bridge\n",
    "from damage import damage\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 输入桥梁类型，格式为'HWB'+'1-10'\n",
    "bridgeType = input('Class of bridge?(from HWB1 to HWB10)') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 获取桥梁易损性函数中位值\n",
    "slightMedian = bridge[bridgeType]['slight']\n",
    "moderateMedian = bridge[bridgeType]['moderate']\n",
    "extensiveMedian = bridge[bridgeType]['extensive']\n",
    "completeMedian = bridge[bridgeType]['complete']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 输入灾害名称，wenchuan,tangshan,taiwan\n",
    "damageName = input('Name of damage?')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 获取灾害强度相关信息\n",
    "PGA = damage[damageName]['PGA']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "C:\\Users\\29968\\AppData\\Local\\Temp\\ipykernel_15908\\2635857658.py:17: UserWarning: linestyle is redundantly defined by the 'linestyle' keyword argument and the fmt string \"r-\" (-> linestyle='-'). The keyword argument will take precedence.\n",
      "  plt.plot(x1,y1,'r-',linewidth=1,linestyle='-',label='slight')\n",
      "C:\\Users\\29968\\AppData\\Local\\Temp\\ipykernel_15908\\2635857658.py:26: UserWarning: linestyle is redundantly defined by the 'linestyle' keyword argument and the fmt string \"g-\" (-> linestyle='-'). The keyword argument will take precedence.\n",
      "  plt.plot(x2,y2,'g-',linewidth=1,linestyle=':',label='moderate')\n",
      "C:\\Users\\29968\\AppData\\Local\\Temp\\ipykernel_15908\\2635857658.py:35: UserWarning: linestyle is redundantly defined by the 'linestyle' keyword argument and the fmt string \"b-\" (-> linestyle='-'). The keyword argument will take precedence.\n",
      "  plt.plot(x3,y3,'b-',linewidth=1,linestyle='--',label='extensive')\n",
      "C:\\Users\\29968\\AppData\\Local\\Temp\\ipykernel_15908\\2635857658.py:44: UserWarning: linestyle is redundantly defined by the 'linestyle' keyword argument and the fmt string \"m-\" (-> linestyle='-'). The keyword argument will take precedence.\n",
      "  plt.plot(x4,y4,'m-',linewidth=1,linestyle='-.',label='complete')\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEWCAYAAACaBstRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAAsTAAALEwEAmpwYAABcI0lEQVR4nO2dd3hURffHP5NGEggJofcioQUIvSgKiog0saCIwCugYu8FyyuoYPuJXUB4bSigCCqiqChFpQqhhd5DL0kIIT3Z3fP7YzYYIGWTbEsyn+e5z+7dO3fm7ObmnjszZ75HiQgGg8FgMDiCj6cNMBgMBkPpwTgNg8FgMDiMcRoGg8FgcBjjNAwGg8HgMMZpGAwGg8FhjNMwGAwGg8MYp2Eo8yiltiuletnfv6SUmmV/30AplaKU8nWDDTWVUn8rpZKVUm87qc7hSqnfCzj+p1Lqbme0ZTDkYJyGwWMopWKVUun2G3fOVsfZ7YhIpIj8mcfnh0WkkohY7fa48iY7FogHKovIkxcfVEp9oZTKsv8GyUqpDUqpngVVKCKzReQ6F9lrMOSJcRoGTzPIfuPO2Y7nPqiU8vOUYU6mIbBDCl5N+38iUgkIBaYB3+fXCypDv4uhlGGchsHrUEqJUupBpdReYK/9s/eVUkeUUufsT+FX5iofpJSaqZRKVErtVEo9o5Q6mut4rFLq2jzaaWRvy08p9SpwJfCR/Wn/I6XUlIuHkpRSPymlHsvH7suVUuuVUkn218vtn38B3Ak8Y6/7EltyIyI2YA4QDtS01zFKKbVKKfWuUuoM8JL9s5W52u+jlNplb/8jQOU65quUelspFa+UOqiUeijnu9uPhyqlPlVKnVBKHVNKTXLHsJ2h9GGchsFbuRHoCrSy768H2qFvpHOAeUqpQPuxCUAjoAnQBxhR1MZE5AVgBfCQvcfzEDATGKaU8gFQSlUDegNfX3y+UiocWAR8AFQF3gEWKaWqisgoYDb2noSILCnIFvvN+j/AQeBUrkNdgQNADeDVi86pBnwH/BeoBuwHrshV5B6gH/o37ID+fXMzE7AATYH2wHWAmQ8xXIJxGgZPs0Apdda+Lcj1+esickZE0gFEZJaIJIiIRUTeBioAze1lbwNeE5FEETmKvnGXGBFZByShHQXA7cCfInIqj+IDgL0i8pXdxq+BXcCgIjT5lFLqLJAKvAe8mDPfYue4iHxorz/9onP7o4e/5otItv38k7mO3wa8LyJHRSQReCPngFKqJtqhPCYiqSJyGnjX/n0NhgswTsPgaW4UkTD7dmOuz4/kLqSUetI+9JRkv7GGop+oAepcVP6Cc0vITP7tuYwAvsqnXB3g0EWfHQLqFqGtySISBgQBnYC3lFL9ch0v6Htd8BvY506O5Hf8ovcNAX/gRI4DB6ajezQGwwUYp2HwVs5PGNvnL8ahn5ar2G+sSfw7Zn8CqJfr3PolbTMXs4DBSqkooCWwIJ9zj6NvvrlpABwrshGabcAqdA+mIPtyOEGu762UUlz4OxT0Gx0BMoFquRx4ZRGJLKrthrKPcRqG0kAIerw9DvBTSo0HKuc6/i3wnFKqilKqLvBQMds5hZ4XOY99uGs9uofxXR7DQjn8AjRTSt1hn1gfip6P+bk4hiilWgA9gO0OnrIIiFRK3Wyf3H4EqJXr+LfAo0qpukqpMLQTBkBETgC/A28rpSorpXyUUpcVFvJrKJ8Yp2EoDSwGfgX2oId8MrhweOUV4Ch64ngJMB/95FxU3geG2KOwcs+LzATakP/QFCKSAAwEngQSgGeAgSISX4T2c6KrUtE38c/Rw0SFYm/nVvRcRQIQge6p5PA/e50xwCa0k7MAOXMm/wECgB1AIvo3rF0E2w3lBGWSMBnKGkqp+4HbRcQpT8pKqavQw1SN7OGwpR77XMnHInLxkJrBUCCmp2Eo9SilaiulrrAPqzRHP+3/4KS6/YFHgU9Ks8Owr2Xpbx86q4sOU3bKb2QoX7jMaSilPlNKnVZKbcvnuFJKfaCU2qeUilFKdXCVLYYyTwB6GCcZWAb8CEwtaaVKqZbAWfQwzXslrc/DKOBl9NDTJmAnMN6jFhlKJS4bnrJ36VOAL0WkdR7H+wMPo+PLu6JjyLu6xBiDwWAwOAWX9TRE5G/gTAFFBqMdiojIWiBMKWUm3gwGg8GL8aToWV0ujIA5av/sxMUFlVJj0SqhVKxYsWOLFi3cYqDBYDC4nZzRH5G8t7yO5Vc+V30bDhyIF5HqJTXPk05D5fFZnmNlIjIDmAHQqVMniY6OdqVdBoOhrGK1QkYGpKfnv2VkQFYWZGbq15ytpPtZWWCxXLplZ//73moFX1/w9wc/v6JtF5/j6/vvq68v6sCBixULioUnncZRLlyVWg+9qtZgMJRXRPSNOyXl0i05Oe/Pc46lpRXsDNLT9Y05MBCCgvLfAgOhQgUICPj3NWerUEEfr1z5wmN5lbv4vb9/oc4gxZpOoH8QVpuVI+eO0DS8KfvO7EOhuCz8Mr7d/i39mvYjKTOJ3/f/zpj2Y5i/Yz7VgqvRq1Ev7vrxLt69/l32JOzh661f83bft3l+6fN0qtMJvv3WKX8iTzqNhcBDSqlv0BPhSfaVqQaDobQiAufOQXw8nD2rt6Skf98Xtp+crG+ylSrlv4WE/Pu+QYN/3wcHO+YMVF6DHEX9mkKmNZNAv0CSMpJIykyiQWgD/jn6Dw1CGxAeFM6MDTN4uOvDrDi0gsNJhxneajhPLH6CO6PupFalWtw2/xb+GvUXH/zzASdTTvJa79e4/st+TB0wlYr+FXl88eP8NOwnVh5eib+PP5eFX8bGExvp2VAvP8q2ZgNQo2INQiuEAnBD8xsI8A2gSZUm/CfqPwA82PlBKgZULPF3zsGV0VNfA73QonKn0HHh/gAi8rFdG+cj4HogDRgtIoWOO5nhKYPBjVitcPo0nDypX+PjC94SEvQNumpVqFIFQkMhLOzf7eL93J+FhuoneF/3pPGwiQ2FIt2SzrFzx4ioGkH08Wgq+lekZfWWTFg+gXE9xrEzbiff7fyO13q/xuO/PU73+t25tdWtVHq9EmfHnWXJgSUsObCEt/u+zWsrXqN34950qN2BF5e/yBvXvsGOuB3EpcbRs1FPoo9Hc1mVywipEMLBxINEVI0g05KJUooA3wCXfl+l1AYR6VTiekrbinDjNAwGJ2C1akdw+DAcPw4nTujt5MkL38fHQ3g41KoFNWtCtWoFb1Wr6qd5D5GUkUS6JZ1alWqxaM8iutTtgiB88M8HTLpmEtOjp5Nty+ahLg/RckpLfrz9R7KsWby64lW+vuVrvtryFdUrVuf6ptczZd0U7mx3JxmWDGLPxtKpTifi0+IJ9g8m2D/YY9+xuBinYTAY8ic1VTuEw4fh0KF/3+dsx45pZ1C/PtStC7Vr661WrQvf16ihx+E9iIiwO2E3Laq1YG/CXvad2Ue/iH688tcrXHfZdbSq3oq209oS+1gss2Jmsf/Mfib0msCkvycxNHIotSrVYuHuhQxvO5yTKSdRKGpWqolNbPio8iOKYZyGwVDeSU6G/fth3z7Yu1e/5mxnzujx/pytYcML9+vV02P8Hibnxv3bvt/o1agXh84eYsGuBYzrMY5nlzxL6xqtGd5mOO2mt2PtXWvZFb+LTSc3Mab9GNYdW0eD0AbUrFiTM+lnqBpc1dNfx6txltMwyekNBm/nzBnYvv3CbedOPYF82WUQEQFNm0LXrjB8uN6vWxd8PP8UneMUlh5YSlStKCw2C6+teI0P+n3ApL8nISK82PNFftj5A+1qtSMsMIxW1XWG36cuf4pg/2CUUmy5bwsA7Wu3p33t9gB0qdvlfDvGYVyI1QpxcVC9uh5pXLrUeXUbp2EweAtZWbBjB2zcCFu2/OsgUlOhVSuIjITWrWHQIL1fp45XOAaAg4kHqRRQibDAMJ78/Unev/59vor5ipWHVzJj0Az+OvQXdULq0DCsIX2a9AFg3BXj8PPRt6Dpg/5VgB/UXGfIrRZc7dKGDNhs+s++ahUcPKhHGn184Omn4dVXYepUHbNQpQrExOipqeXLnde+GZ4yGDxBVpb+j964ETZs0K/bt0PjxtChA0RFaScRGannHZwQJlpSRASlFLNjZnNTy5vYdnobU9dP5Ysbv+DFZS/SvX53+kf0Z8aGGdzV/i6sYsXPx69czRuUFIsFjh6F2Fg9FXX11XqJx+jR+vNjx2DUKHjnHbjzTr0usG5daNkSxozRx202PR118VSUmdMwGEoTp0/DmjWwerXeNm2CJk2gUyftJHIcRUXnxdMXF4vNgoiQlp3GJxs/4cnLn+SdNe+QmJ7IxGsm8vzS53m066NUCqhEXFocjcIaedrkUkN2tl6OUr06/PGH7i3kOIjFi2HOHBg/Xk9BNWoEzzyjRxuXLtXOoV49HaBWnGcIM6dhMHgzR47AsmV6XGDVKh262q0bXH45TJgAXbroNQkeJsuahb+PP7vidxF9PJqRUSO5dd6t3NvxXno16oXFZgHgvk73UcFXh9K+1vu18+c7c9FYWcFi0U5g7179xN+uHQwdqjuThw/DddfBTz/pnoPNBlddpR2Ej4/uRYwadWmdAwZc+pmnMD0Ng8EZxMXBn3/qR8JlyyAxEa65Rm89eujxAy+Zf1h1eBWC0KNBDxq/35g1d60hOTOZ6OPRDGszDIvNcn6uwZA/iYl6ZHHXLu0gxo7VAWmRkdpZRETAXXfB7bfDX3/pZS6NG3tuGYvpaRgMnsRm04+OP/0EP/+sw1yvvBJ694YHHtAT1l7gJKKPR1MnpA4VfCvQa2YvYu6LId2Sjs2ehHDPQ3vw9/WnVqVaRFSNADAOIw/WrYOVK3Wcws6d/z4bfPghtGihHURIiB4+Onv20mjmnk5JPOwdmJ6GweAoaWmwZIl2FIsW6eGlQYNg4EA97OThRXApWSkE+QWxPW47s2Jm8X99/o83Vr5BjwY9uKL+FRxPPk7dynU9aqM3kyMw+7//6emnrVv1zX/FCh2RtHu3Dlpr2RK6d/f4n7vImIlwg8EdpKZqJzF3rn687NTpX0cREeFp69iTsIcTySfo2agn3T/tzqc3fEqdkDrsiNvB5fUv97R5XonN9q8i+Rtv6OjmzZv1n3TqVHjrLS2D1bq1dhJhYZ622DkYp2EwuIqMDN2TmDtXh7RcfjncdhvceKMOfvcgFpsFH+XD7fNvZ+aNM9l4YiN7EvYwuv3ocieL4Qg5t7eDB2HKFD3MtHmzdg7DhsHEidC2rQ5ca9LEK0YUXYZxGgaDMxGB6Gj4/HPtLNq10zOYN9+sYxw9hMVmYeOJjXSp24W3V79NSlYKE3pN4I/9f3BlwysJ9PO8FIi3ERMD332nHcT69XpEMTQUvvlGB6116OBx3+8RjNMwGJxBXBx8+aV2FunpOt7xzju1PpMHWXFoBSEVQmhetTmDvh7EbyN+I8uaRaBfoOlN5GLnTu0UVq3SS1+2btVBbH/9pR1E58564Xx5J+tUFhVqVTDRUwZDsVm/Hj76CBYuhBtu0GMXV17p0fGJM+lnmPT3JN7p+w6nU09jsVkI8g9iyX+WACaqKS3t3yim+Hh47z1YsAAOHID+/fVQk7+/Xgdx3XWettYzWFOtpGxNIXVLKilbUsiOyyZyXiSn5552Whvl+yo0lC8yM3XKy48+0iu0H3hA6zF4aPgpR5bj/1b9H43CGjGk1RCiakZhExu3tLrFIzZ5E+npugeRmakXt/Xvr9VXevSAPlq+iuee86yNnkREOLf6HJUvr0zCogQOPHOAjNgMglsGUymqkt6GVgKg3iP14FHntGuGpwxln+RkmDED3n1Xx0s+8oi+A7kpQ1xuRASLzcLZjLNcN+s6NozdwIHEA1QNqkqVoHI40J6L7Gz9p6pSBfr10z2KqCgYMQLuv/9fob7yiOWchZRNKSRvSCZ5QzLh14VT8z81iekXQ+S8SKznrGTHZxPcKhgf/7x/JDOnYTAURlycXn01bZpedDduHLRv7xFT0rLTCPYP5pk/nqFpeFPGdhzLkaQj1A+t7xF7vIX0dL0uYvFi7SSeeEKrrKxZo0NeQ0I8baH7EatgSbLgH+7PgecPEP9jPBmxGVRsU5GQjiGEdAohrFcYQY2DilSvWRFuMORHfDy8+SZ8+incequ+AzVt6nYzbGIjw5JBtjWbDjM6sOvBXUzoOeF8qtDy6DAyMuDvv+GXX/Qq6YED9SrrUaNg5kydMRb04rnyQnZCNufWnuPc2nPUf6Y+iX8kcubXMzT/X3NCrwql+pDqVGxTMd8ehLsxTsNQdjh3Tg9BffihXlexbZtHQmdSs1KpGFCRN1a+ga/yZVyPcex4YAf+vv74+5ayZcRO4PBhCA7W6yTbtNHrIvr31+/9/eHjjz1toXuxJFuI+y6OpL+TSFqZRNapLEI6hxDaPRTJEqrdVI3qN1cHoOr13pdcyjgNQ+knI0Ov1nrzTejbV4fYNGnidjMsNgs2sdH247ZsuncTT3Z/kgDfAAAq+HlIpc5DxMXBJ5/A99/rhXVffQXXX69lwMPDPW2d+xARsk9nE1AzgEOvHSL0ilCCI4NJXJxI6FWh1Hu8HhVbVUT5ej5fiqMYp2EovYjAjz/qgfA2bbTMR+vWbjZBR0DN2z6PpQeX8vHAj9nxwI5y5yRE9DqJH37QI4JhYToh0JtvaulvP/udpqw7DLEJqdtSSVyWqHsSK5Ko1K4SUX9EEd4vnAp1KhBQLYBWX7fytKnFxkyEG0onO3fCo4/qO9P778O113rEjIFzBvJ679e5LPwyfJRPuVqhLaK3zZv1wnk/P7jlFrjvPi0BXl7IPJ5JhToVOD3/NHsf3ItfZT/CrgkjrGcYoVeGEljfO64JEz1lKJ8kJcHLL+vxjhdf1LGYbpYb3X9mP5NXT2bawGkcTjpM/cr1UV6QjtVd7N4NX3+ts8xNmaJzS8XG6k5eefgZss9mc/bPs/iH+xPaI5T1rdfTfnV7bOk2bJk2ghoVLarJXZjoKUP54+eftZPo21eH3FSv7tbmF+5eSJXAKnSr141bI28FoEGoZ+VG3MXx41q/aeNGnYXu9tth9mwt+quUHh0sq4hNSN6QTMKiBBIXJ5K6LZXKl1emzv11UD6Kzts764eGME9b6h6M0zB4P/Hx8NhjsHat1om6+mq3NS0i/L7/d6677DqC/YMJ8g/C39efaxpf4zYbPEVmplaF/+wz/dMvXKgFf48c8ci6SLdiTbOSEpNCaLdQdt+1m3NrzxE+IJzGkxpT+YrK+Ab++wOUp14mGKdh8GZEYP58vYJ72DCd+KCi+3JSJ6YnEhoYyuyts2lfuz3XNvHMvIm72bJFRyofOaKD0saM0X+G4GBPW+ZaUnelkvBzAlWurUJAzQCOvX+M0G6hREyLuMBJlHfMnIbBO0lKgnvv1TrXn32mB87dyLpj63j5r5dZdMcit7brKVJTYdYsvTr79Gn9/qqrPG2VaxERktcnE/9DPHE/xGFNsVJtUDXqPlSXipHuezhxF2ZOw1B2WbsW7rhDrwDbuPHShMsu5MN/PiQ0MJSRbUeyYOgCt7XrKXbt+rcH8dtv8OqrOhCtrA4/2Sw2sk5kEVg/kC29t5B5PJPqN1Wn5VctCekUUu6GmoqDcRoG78Fm04H9770H06frTHluQET4csuX3BZ5G/0i+lElsApKqTK7etti0ctbpk6F7dt172LQIL3Goixiy7aRuCSRyl0rkxKTwqmZp2jxeQsi50fiH142/8auxDgNg3dw4gSMHKm1r6Ojob57dJnOZpwlLDCMQ0mHSEhPoGm4+zWq3EVioh71q1JF++R77tHrKwICPG2Z8xGrcPavs5z+5jTxP8QTFBFEsxnNqNKrClV6aTVh4zCKh3coYBnKN2vW6BRrV1wBy5a5zWGcSD5Bry96YbVZGd9zPPUq13NLu+4mNlYHn112mY6ACg2F33/XYbNl0WHsf2Y/a+qtYf9T+wlqGkTH6I50WN2BSq0redq0MoHpaRg8y4wZ8N//6snugQPd0uS6Y+tYc2QNj3Z7lHX3rMPXp2wO4Kek6PmKG27Q+SliYqBe2fSLnF15ljO/nKHJa02o1KESte+qTXDzMh7u5SGM0zB4huxsePhhrZO9ciU0a+byJuNS41BKUa9yPVrX0BpVOYKCZQUR3VmbOFEnJPzuOy3zUdaSF2UnZnP6m9OcnHmSpu81pWLLiviF6ttZzdtreti6so1xGgb3k5QEQ4ZAhQrwzz9uy7QzZf0UmldtzrA2w6gT4n7JdHcwerQOPnvhBb20BcqOw7BZbCQuTuTkzJOc+f0M4X3DaTS+ESGdQvDx88G/qpmjcAcuXaehlLoeeB/wBT4RkTcuOh4KzAIaoB3YZBH5vKA6zTqNUk5srE74fM01OveFn+ufW8b9MY6RUSOJrB5Z5kIqRWDRIi3FNXu2lvuoW7dshcxa0634BvmyscdGsELNO2tSY2gN/KsYJ1EUvH6dhlLKF5gC9AGOAuuVUgtFZEeuYg8CO0RkkFKqOrBbKTVbRLJcZZfBg2zcqGM7x43Tq7xdzIHEAzSp0oT+Ef1pENqgzDmMv//WE9xWq9Zu9PGBBmVECsuWZdMy41tTOfj8QaL+iKLtorbnh6AMnsOVHdcuwD4ROWB3At8Agy8qI0CI0v/NlYAzgMWFNhk8xV9/6Sw8H37oFodhsVkY8f0ITqWcomejnlSuUNnlbbqLVat0bm0/Pxg/XuexGDKkbAxDpR9M58BzB1jTYA1nfjlDSMcQ2vys1RCNw/AOXHmZ1QWO5No/av8sNx8BLYHjwFbgURGxXVyRUmqsUipaKRUdFxfnKnsNrmLhQp2Z55tv9MIAF5KcmcwLS18AYNWYVdSsVHYmRWNidIDZ8OGwf78WD7zxxtLvLGwWG3EL4ojpF8OGzhuwZdpo92c7qt9cHeWj8KlQyr9gGcOVrjuvsYCLJ1D6ApuBa4DLgD+UUitE5NwFJ4nMAGaAntNwvqkGlzFrFjz1lB5479zZpU1lWjIJ9g+mYVhDoGypjx48CNddB889pyOiKpSBxIC2bBs+/j7suXcPaTvTqHN/HSK/j8Q3qAxNyJRBXOk0jgK5V2nVQ/cocjMaeEP0bPw+pdRBoAWwzoV2GdzFzJnw/PM6BrSVa9Nbxp6N5bZ5t/HP3f8wtuNYl7blLs6dg0mToHJlvZQlNtatMlwuQ2yCNdVKdPtoumzvQsSHEfgGG0dRWnBlv289EKGUaqyUCgBuBxZeVOYw0BtAKVUTaA4ccKFNBnfxxRc67nPpUpc6jPTsdP6K/YtGYY34feTvZaZ3MXMmtGihU4ncfbf+rDQ7DBHhzJIzxAyMYf+T+/EL8aPThk74VPAxDqOU4bKehohYlFIPAYvRIbefich2pdR99uMfAxOBL5RSW9HDWeNEJN5VNhncxOef6xnapUuheXOXNnU46TDf7fyOno16EhYY5tK23MGRI1pFJSFBiwq6eETP5diybJyafYqj7xxFRKj3WD1qDtfzTGZiu3Ri8mkYnMvcufDEE7B8uUtXeW8/vZ15O+bxUq+XXNaGOzl6VEcir1mjM9mW5l5FDmdXnmXnsJ0Etwym/tP1qXJtlTLTEyyNOGudhglLMDiPRYvg0Ud1YgYXOoxsazb1Kteja92uLmvDnSxdCu3aQZMmOkKqtDuMox8eJWVrCsHNgmm9oDVRv0cR3ifcOIwygnEaBufw119aw+LHH6FNG5c1s/TAUv6z4D+EBobSL6Kfy9pxBxs2wM6d0LEjrF6t9aIqlVIh1sxjmRx5T0fYB0UE4VfZj4AaAYR0dI9EjMF9GKdhKDmbNv27DqOra57+rTYrBxMPcnXjq5nSf4pL2nAXKSnw+OM6MeHBgxAW5ha9RpeQeTyTvY/sZX2b9WQezsSWbaPq9VUJbFjKu0uGfDFOw1AyDh/W0iBTp2o9KRex7tg6Xlz+Ij7Kh/CgcJe142pEdDrVhATYtk07jtJI5vFM9j66l/Wt16MCFF12dqHpO03x8Te3lLKOCV8wFJ+zZ/Vd74kntI6FC0jKSGLl4ZUMaDaAbvW6uaQNdxAfD1Om6PUWixZp2fLSiNgEFGwdtJWwXmF03tGZCrXKwEpDg8OYxwJD8cjKgltugauv1mMtLiIuLY5/jv0DlN4V3j/+CG3b6sV6VmvpdRipO1LZ0mcLSik6/NOBpm83NQ6jHGJCbg1FR0SvOIuPh++/d4kO9/Hk40yPns5LvV4qtc4C9AT3f/6j1zr26OFpa4qONc3K0Q+OUqF2BWr+pyZZp7KMoyiluE0aXSnVWkS2lbQhQxnio49g/Xp9R3RR4obKFSrTrGqzUuswFi/W8xbDhum5i9IWRmvLtnHys5PEvhJL6BWhVJ9UHaWUcRgGh4anPlZKrVNKPaCUCnO1QQYvZ9kyePVVPebigvjQhLQEhnw7hADfAIa3He70+l1NRobOcXHPPVC7NihVuhyGiHB6/mnWR67n9LzTtF7QmshvIwluZvJtGzSF9jREpIdSKgIYA0QrpdYBn4vIHy63zuBdHDwId9wBc+ZA48ZOr94mNsKDwnm066OlNnf3U0/ByZM6L3d4KQvyEhGSo5M5NOkQEVMiCO9Tyr6AwS04PKdhz8R3I/ABcA6tFfW8iHzvMuvywMxpeIi0NOjeXT9CP/SQ06vPtGRy5edXsnjEYqoEVXF6/a5EBD75RGexDQ2F4GDdwyhNHJt6DFu6jfpP1kdsgvIpZV/AUCjunNNoi5YwHwD8AQwSkY1KqTrAGsCtTsPgIR56SIcAPfig06vOsmZRwa8C3976balzGGfP6oXwhw/rZSp16njaIsexplk58s4Rao+pTfVbquMTrEerjcMwFIQjcxofARuBKBF5UEQ2AojIceC/rjTO4CV8/jmsXQvTprnkEXrAnAFsObmFRmGNnF63K7Fa4aqrtCrt6tVw2WWetsgxRITTc0+zrsU6UremgkBAzQD8QsyyLUPhOHKVfC8iX+X+QCn1qIi8f/HnhjJITAw884zWlnLyxHdqVipB/kHMHTK3VK3yFtEig71764V69esXfo63kLozlb0P7SU7PpuWc1oS1iPM0yYZShmO9DT+k8dno5xsh8EbSU7WmlLvvuuSREovLn+R2TGzS5XDSEmBESP0IvikpNLlMOK+i2PzVZupNrgaHTd0NA7DUCzynQhXSg0D7gB6ACtyHQoBrCJyrevNuxQzEe5GRo/W6zA++cSp1WZbs0nNTqWCbwUq+FXAR5UOYYJTp/S8RbdueqlKUJCnLXKM+IXxBLcIxi/cD8kWKtQ2ay3KI+6YCF8NnACqAW/n+jwZiClpwwYvZ948WLUKNm50etXfbv+WDSc28E7fd5xet6s4dw6qVYM33tD6jKWBrFNZBNQMwJpsxZpiNWstDE7ByIgYLuXoUejQAX7+Gbp0cWrV8WnxVA2qisVmwd/X36l1uwIReOstmD1bK8D7lIJOkc1i49j7xzj0+iE6rutIUJNS0iUyuBSXZ+5TSq20vyYrpc7l2pKVUudK2rDBS7HZtFjSI4843WEcOnuIAXMGIEipcBhpaXot47x52n+WBoeRsiWFTd03kfBLAh3/MQ7D4HzyHZ4SkR72V5N6qzzxzjuQnQ3PPefUahPSEmgY1pCVo1eWmjmMI0f0Yr2///b++QtrhpVDEw9x4n8naPJGE2qNrlVqdbsM3k1BPY3wgjZ3GmlwE7t26UH7mTOdKkRoExv95/TncNLhUtHD2LxZy4E0bw4ff+z9DgMgfkE86XvT6RTTidpjahuHYXAZBU2EbwAELRdyMQI0cYlFBs9gtepoqZdegibO+9MmZSRRKaASK0avKBV6UosW6Z9hSinJKHvk3SME1AygxrAa1Ly9pqfNMZQDChqecr4incF7ef99CAiABx5warUT/55IVM0oRkaNdGq9rmDNGi2ttXChDqv1Zs6tP0dwy2CqDqyKf1V/07MwuI2C1mm0EJFdSqkOeR3PkRNxNyZ6ygXs2QOXX66lQpo2dVq16dnp+Pn44evj69XzGDYb7N0LERFw+jTUquVpi/LHlmUj9uVYTnxygjY/t6Fy58qeNslQSnDHOo0ngLFcuEYjBwGuKWnjBi/AZoO77oIXX3Sqw1h/bD3j/xzPr8N/dVqdriAzUweLZWToFCHe7DBSYlLYOXIngQ0D6bSlk0mIZPAIBQ1PjbW/Xu0+cwxu59NPdbSUE+XOrTYrnet2Zu6QuU6r0xWcOwc33QRhYTDXi00Vm3Bk8hGOvHWEJm81odadJjLK4DkKHTNQSgUqpZ5QSn2vlPpOKfWYUqoU5SIz5Mvp0/DCCzB9utOipaw2K90+7cbx5ONUruDdQyenT0OnTvDtt96dXS8jNoPEpYl0jO5I7VEmMsrgWQpdEa6U+hYtHTLL/tEwoIqI3Opi2/LEzGk4kZEjoWZNmDzZKdWJCEopTqacpFYl7x3n2bdPa0e9+653J0tK3pjMqVmnaPqO84YNDeUXtyVhApqLSFSu/eVKqS0lbdjgYZYu1avWtm93WpWvrXiNOiF1GN1+tNPqdDbbtkHfvjB+vPc6DGu6laxTWQQ3D6bG0BqeNsdguABHnMYmpVQ3EVkLoJTqCqxyrVkGl5KRAfffDx9+6NQcGfd3vh9v1jI7ehSuvVb3MIYN87Q1eZOyLYWdw3ZSdVBVmrzWhMpdvXuIz1D+yNdpKKW2oqOk/IH/KKUO2/cbAjvcY57BJbz9NrRsCTfc4JTqTqWc4snfn+TLm7702tDatDSoWxd+/11nrfU2RIQTn57g4HMHz092GwzeSEE9jYFus8LgPo4e1fpS69c7rcqqwVW5q/1dXuswli7Vi/a2bvVOh2FNtbLngT0kb0im3Yp2VGxR0dMmGQz5UlDI7aHc+0qpGoAXx5gYHGLcOLjvPqdJhXwc/TFX1L+Cqxt7Z2R2jizI/PlQ0UvvxbtG78In0IeO/3TEt6LzNL8MBldQ6JyGUuoG9AK/OsBp9PDUTiDStaYZnM6qVTrX9/TpTquyRsUaXpuuNTMT/vtf+Okn6NrV09ZcSvzCeMKuCaPZ9Gb4hfmZUFpDqcCR8YSJQDdgj12PqjcOToQrpa5XSu1WSu1TSj2bT5leSqnNSqntSqm/HLbcUDSsVp0j4803nTL5nZCWwNxtc7m55c3UrVzXCQY6l3XrdP6L6GjvdBgAyeuTyTqRhX8Vox1lKD044jSyRSQB8FFK+YjIcqBdYScppXyBKUA/oBUwTCnV6qIyYcBU4AYRiQQ8svajXPDFF3oF2x13OKW6+LR4jicfd0pdzua332DgQC2p5USFd6eQeSyTzddsJv1AOo0nNiY4wqRgNZQuHAm5PauUqgSsAGYrpU4DFgfO6wLsE5EDAEqpb4DBXBh5dQfwvYgcBhCR00Ux3uAgKSlaW2rhQqcsTth0YhORNSJ5vPvjTjDOuSxerLWkfvwRIr1sADVpVRLbb9tO3YfqEtjYTA8aSieO9DQGA+nAY8BvwH5gkAPn1QWO5No/av8sN82AKkqpP5VSG5RS/8mrIqXUWKVUtFIqOi4uzoGmDRfwzjvQq5fWzHACH637iN3xu51Sl7NZvx5++AG6d/e0JRdyfMZxtt20jeb/a07D5xqa4ShDqaXQnoaIpCqlaqF7DmeAxfbhqsLIL3nTxe13RM+TBAFrlFJrRWTPRTbMAGaAlhFxoG1DDqdO6VwZTpBeERESMxL5dPCnTjDMuaxcCRaLnvj2JmxZNvY+spekv5Nov7I9wc3McJShdOOIYOHdwDrgZmAIsFYpNcaBuo8C9XPt1wMuHgQ/CvwmIqkiEg/8DURhcB4vvwx33gmNS55Ta+3RtYz+0fskQjZuhJtv1mK93oZkC77BvnRY28E4DEOZwBHBwt3A5Tm9C6VUVWC1iDQv5Dw/YA+6F3EMWA/cISLbc5VpCXwE9AUC0M7pdhHZll+9RrCwCOzeDT166NzfVauWqCqb2PBRPmRbs70qz/fOnXDNNTB1qpY59xay4rPYOWInbX5qg4+/dy56NJQvnCVY6MjVfBStcptDMhfOVeSJiFiAh4DF6HUd34rIdqXUfUqp++xldqLnSWLQDuOTghyGoYg89xw8/XSJHQbAzXNvZsPxDV7lMAD8/OCDD7zLYWSeyCSgWgBNXm9iHIahzFFQutcn7G/bAW2AH9FzEoOBdSJynzsMvBjT03CQdevgllt03GlQUImrO3buGLVDanuNVMjx4/DGG/Dee3o9hjcgIhz5vyOc+OQEnbd3xifASwwzGHBPTyPEvu0HFvDvJPaPwImSNmxwMf/9r95K6DD2ndnH4789Tt3Kdb3GYZw9q+XNa9b0Hodhy7Kx++7dnJ57mnZ/tjMOw1BmKUh76uXc+0qpEP2xpLjcKkPJ+Osv2L8fxjgSr1AwNSvW5IbmzlHDdQZZWTB4sJ7HeP55T1ujyU7MZvst2/EN8aXd3+3wq+TI8ieDoXTiSPRUa6XUJmAbsN2+nsLLlk0ZziOiF/JNmAD+JZt/+G3fb5xKPeVVYoT+/vDww96TdS/jSAYbu2+kUvtKtP6+tXEYhjKPI33oGcATItJQRBoCTwL/c61ZhmLzxx8QFwfDh5e4qiNJR0jOTC68oBsQ0QK9a9bAkCHeMSwlVsG/mj+NJzam6dtNUb5e4MUMBhfjyL9eRbveFAAi8ifgpSLT5RwRPY/x8sslFl3afno793S8h/a12zvJuJLxf/+nNaW8RRok43AGm67chE8FH2rcalKyGsoPjjiNA0qpF5VSjezbf4GDrjbMUAx+/lnrgQ8ZUqJq4tPieejXh8i2esdqua+/hmnT4NdfITTU09ZA2r40AhsEEvldJMrH9C4M5QtHnMYYoDrwvX2rBnjfsuDyjgi88gqMH1+isZssaxZVg6qy/M7lXrMmo3t33cuoU8ezdogIB188yLYbtmHLtlGhdgXPGmQweIACZ+3s8ubzRORaN9ljKC6LF0N6eolXub239j0AnrniGScYVTJ274bXX4fPP/f8pLfNYmPPvXtIjUnVIbVm0Z6hnFKg0xARq1IqTSkVKiJJ7jLKUEREYOJEeOGFEs8QP9n9SdIt6U4yrPjEx8OAAforedphWFOtbB+6HawQtTzKREgZyjWO3GEygK1KqU+VUh/kbK42zFAEli/Xd9nbbit2FSLC3Qvv5njycSoFlDyzX0nIzoYbb4Rbb9X5vT1qS2I2m3tvxr+aP60XmpBag8GR/4BF9s3grUycqFe6lTBiamjkUOqEeHjiAK0n9eyz0L+/py0B32Bfat9dm9p31TY5MAwGCp/TuBE9Cb5VRBa7xSJD0Vi5Eg4dKlEa13OZ5/h5z8/c0cY5qWBLwhtvQMeOOl2rJ7GmWdk6aCutF7Smzt2ed6QGg7eQ7/CUUmoq8DhQFZiolHrRbVYZHOfVV/VjeQlWf8enxXM46bATjSoe332nJc7btPGsHZYUC77BvjR9pyl+IWY4ymDITUEqt9uAKPtkeDCwQkQ6utW6PDAqt7mIiYF+/eDAAahQvPDPkyknqRJYhQp+ng0f3bpV60n99pvuaXiKc9Hn2HbjNjr+05EKdU1IraHs4A6V2ywRsQKISBp5p281eJL/+z949NFiOwyA6dHTmRUzy4lGFY/wcJg507MO4+xfZ9nafyvNpjYzDsNgyIeCehppwL6cXeAy+75Cq922dYuFF2F6GnYOHYIOHXQvo5jLpEUEpdT5V09gsej1iM8+C5Ure8QEABJ+SWDXqF20+roVVXpX8ZwhBoOLcFZPo6AB25YlrdzgQt55B+6+u9gOw2qzcvlnl/PzsJ+pXrG6k41znKef1ilbK3pQzSxuQRx77t1D64WtCe3mBTolBoMXU1A+jUPuNMRQBBIS4KuvYFvxM+P6+vgyd8hcjzqML77Qclnr1pU4WrjYiAgJPyXQ9pe2hHQM8YwRBkMpwmghlEamTIGbby62GNPJlJO8ufJNGoU1cq5dRUQEfvwRqnhoNCju+zgyDmbQ4tMWxmEYDA5i4glLG+np2mn89VexqxAR6lWu50SjisaJE/D3355f7W05a8GaZvWsEQZDKcM4jdLGrFnQpQu0aFGs02PPxlLRvyLD25Y8SVNxyM7Waid9+nikeQBOzjyJT7APtcfU9pwRBkMppVCnoZTaClwcYpUERAOTRCTBFYYZ8kAE3n9fb8Xkj/1/YBMb93a614mGOc5zz0FIiM4V5QmO/+84h145RNSSKM8YYDCUchzpafwKWIE59v3b7a/ngC+AQc43y5AnS5dqyddrrinW6RabhXs63uNkoxwnK0sPTX31lWfStR6ffpxDrx4iankUwU2D3W+AwVAGcMRpXCEiV+Ta36qUWiUiVyilRrjKMEMevPcePPZYsbXCr5l5DTMGzaBFteINbZWEgwd1WO3s2W5vGoATn53g0CTjMAyGkuLI814lpVTXnB2lVBcgRzvb4hKrDJeyZw+sX18iYcLvh35P86rNnWiUY6Slaanzn35ye9MAnFlyhoMvHiRqqXEYBkNJccRp3A18opQ6qJSKBT4B7lFKVQRed6Vxhlx88AHccw8EBRX51JSsFB759RGqBFbxyMrvhx6C1q1hzBi3N401w0rYVWG0X9me4GbGYRgMJaXQ4SkRWQ+0UUqFomVHzuY6/K2rDDPkIjFRj+ts317sKno06IGvj/tX0O3dC1u26Ahhd/ur9Nh0dgzdQYe1HQhqXHRnazAYLsWR6KkKwC1AI8Av50lVRF5xqWWGf/n0U537tBiL+eLT4jmVcorbIouf1a+4JCZCRIRnVnxnnsgkqFEQbRe3NcmTDAYn4sjw1I/AYPT8RWquzeAOLBb48EM9AV4Mtp3expytcwov6GSSkqBTJ884jPif49nQcQPZZ7PxDyt+nhGDwXApjkRP1ROR611uiSFvfvwR6tfXd+AiYrFZ6NWoF70a9XK+XQUgoqdfrr9er0N0J2f+OMPuMbtp83Mb4zAMBhfgSE9jtVLKw7nUyjHvvadzZhSDp39/2iO5Mn74AXbvhrffdm+7SauT2HnHTiK/i6RyFw/qrBsMZZh882mcL6DUDqApcBDIxOTTcB+bNsENN+hFDn5FV3xJz07HJjYqBrhXd9xigbg4qO1GlY7kzcnEXBdDiy9bUPX6qu5r2GAoJbgjc18O/YAI4Dr06u+BmFXg7mHaNLj33mI5jJf+fImE9AS3Ooy0NOjb1/0OA+D4x8eJmBJhHIbB4GLyvRsppSqLyDkg2Y32GHJISoJ583SGoiIiIjSp0oRqwdVcYFj+PPEEVK0KtWq5r82MwxnYMm00m9bMREkZDG6goEfYOehexQa0YGHu/0gBmrjQLsOXX8J11xXrDvzPsX8Y2XakW2+i8+bBH3/oETV33ruTVidhTbISHGEW7hkM7iDf4SkRGWh/bSwiTeyvOZtDDkMpdb1SardSap9S6tkCynVWSlmVUkOK/hXKICJ6aOqBB4p8akJaApP+noRV3Jsnonlz7Tjclec7OzGbuAVx1Ly9JnXuLV4yKoPBUHQKGp7qUNCJIrKxoONKKV9gCtAHOAqsV0otFJEdeZR7E1jsqNFlnpzl01ddVaTTRITwoHB+vuNnFxl2KdnZOkrq8cehQgX3tGlJsbC1/1YqX16Z6jd6Ll2twVAeKWgi/G37NgX4B5gB/M/+/gMH6u4C7BORAyKSBXyDXiR4MQ8D3wGni2B32WbaNLj//iKP86w8vJLbv7u98IJOZPx4nYXP301LIqwZVrbduI3gyGAum3yZexo1GAznybenISJXAyilvgHGishW+35r4CkH6q4LHMm1fxTomruAUqoucBNwDdA5v4qUUmOBsQANGjRwoOlSzIkT8PvvMGNGkU/t0aAHl4W770a6YgXMnAmbN7snP4ZYhZ3Dd+If7k/z6c3NxLfB4AEc+VdvkeMwAERkG9DOgfPy+o++eFHIe8A4kYIH4EVkhoh0EpFO1auX8eGITz/V+VBDQ4t02t+H/ua3fb9RJ8R94/sbNsD06VCjhuvbEhH2PLAHS5KFll+1RPkah2EweAJHFgDsVEp9AsxC3/RHAI7EgR4F6ufarwccv6hMJ+Ab+xNjNaC/UsoiIgscqL/sYbHou3AxEk/4+/hTwc9NkwrAvn3FlsMqFmIVAmoGcNnky/Cp4IG0fwaDAXDMaYwG7gdytCz+BqY5cN56IEIp1Rg4hk4Te0EGIRFpnPNeKfUF8HO5dRgAixZpnal27Yp02uGkw3Sq0wl/X/dMLCxapJVNdu50z1zGyS9PEtI5hMavNC68sMFgcCmO5NPIUEp9DPwiIrsdrVhELEqph9BRUb7AZyKyXSl1n/34x8U1uswydaqeAC8ik1dPpk+TPgxq7vqF+vHxMHYszJnjvslvFaBM78Jg8BIc0Z66AXgLCBCRxkqpdsArInKDG+y7hDKrPbV/P3TrBkeOQGBgkU8XEbdMDD/2mFY1mTzZ5U2R+GciGbEZ1B7lZk0Sg6EM4k7tqQno8NmzACKyGZ2QyeBMPvsMRo4sksMQEfrP7s/hpMNuiyR69VWYNMn17SRvTmbHbTsIbFh0B2owGFyHI07DIiJJLrekPGOxwOefw113Fek0pRRvXPsG9SvXL7xwCTl6FHr3hoCAYnWEikT6gXS2DthKxNQIqlxdxbWNGQyGIuGI09imlLoD8FVKRSilPgRWu9iu8sWvv0KjRhAZ6fApVpuV6dHTaVW9lct7GTYbjB4NV1/t+nmMrNNZxPSNoeF/G1JjiBtieQ0GQ5FwxGk8DESic2nMAZKAx1xoU/njk0/g7ruLdEpyVjJHzh3BV7k+l+q0aXDuHDybr3qYc7CkWIjpH0ONYTWoe39d1zZmMBiKRaET4ecLKlVRRDyeG7zMTYQfP657GEeOQKVKDp1iExvJmcmEBhZtAWBxWbJERwI3b+7admzZNk5/fZqaI2ua1d4Gg5Nx20S4Uupye/a+nfb9KKXU1JI2bLAzcybceqvDDgNg44mN3Db/NhcapbFY4H//08NSrnYYB/57gPR96dT6Ty3jMAwGL8aR4al3gb5AAoCIbAGKJr9qyBubTcuGFHFoqlOdTvw8zPVKtm++Cd9+6578GKE9QqlQz30r2g0GQ/FwKI+oiBy56OnPvckayip//QXBwdA5X63GS1i8bzH7zuzjwS4PutAwnUzp/fe1vpQrxQhPzT6FLcNG7bvMWgyDoTTgiNM4opS6HBClVADwCI5pTxkKI2cCvAiP8q2qt6JqsOvzYP/4o86TUd+F0byJyxLZ9/g+2i1v57pGDAaDU3FkRXg14H3gWrRy7e/AoyKS4HrzLqXMTISfOQNNmsCBAxAe7tApO+J2EFohlLqVXRtZlJCgc327kpRtKWy5Zgut5rYyazEMBjfgtolwEYkXkeEiUlNEaojICE85jDLFrFkwYIDDDgNgxaEVrDqyyoVG6YRKXbvqSXBXkXksk60DttL0vabGYRgMpQxHoqdmKqXCcu1XUUp95lKryjoiOiypCBPgNrFxb6d7uS3SdVFTyckwahS8+67Wl3IFlhQLMQNiqHNfHWreUdM1jRgMBpfhyBRnWxE5m7MjIolAe5dZVB5Yvx7S06FnT4dPuW3ebaw4tMKFRunUrddcA4NcKJbrG+RLw+cb0uDZMp6B0WAoozjyPOmjlKpidxYopcIdPM+QH598onWmihCWNGPQDEICQlxoFDz3nGt1pfaP20/NETWpcZuRBzFcSHZ2NkePHiUjI8PTppR6AgMDqVevHv4u0vxx5Ob/NrBaKTUfnbnvNuA1l1hTHkhJgXnzYPt2h0+ZvHoyd0bd6bIkS/Hx8Mgj8OWXrhuWAgjvF05gY6Naa7iUo0ePEhISQqNGjczizhIgIiQkJHD06FEaN3ZN0jJHkjB9qZTaAFyNjp66WUR2uMSa8sC338JVV0Edx3J5iwi+ypeQCq7pZYjAffdB48aucxgnvjiBLc1G3QeMnpQhbzIyMozDcAJKKapWrUpcXJzL2ij0NqGUuktEPgW22/d9lVITRORll1lVlvnkEz0O5CCHkg7xePfHXWbO7Nmwa5cO5nIFZ34/w4FxB2j3VzvXNGAoMxiH4Rxc/Ts6MqjeWyn1i1KqtlKqNbAWcO3gelll+3aIjYV+/Rwqfi7zHLd8ewtZ1iyXmeTnB1995Zq5jJQtKewcsZPI+ZFUbFHR+Q0YDAa348g6jTuAmcBW4BfgMRF5ytWGlUk+/VQnpnBwHKhyhcpE3xNNgG+A002x2eCXX2DoUGjvgli4jCMZbB24lYiPIgi7Msz5DRgMLqZXr17kLCTu378/Z8+edbh8bjZv3swvv/ziChM9giPrNCKAR4HvgFhgpFIq2MV2lT0yM/Uj/ZgxDhXfk7CHEd+PcFlXc+pUeOUVsLpARcySYmFr/63Ue6yeiZQylAl++eUXwsLCinVuuXMawE/AiyJyL9AT2Ausd6lVZZEff4S2beGyyxwq3iisEU92f9IlpuzeDS+95LpoKd9gXxq/2ph6T9RzfuUGgwtITU1lwIABREVF0bp1a+bOnXvB8UaNGhEfHw/AxIkTadGiBX369GHYsGFMnjz5fLl58+bRpUsXmjVrxooVK8jKymL8+PHMnTuXdu3aXVJvacSRW0YXETkHIFqo6m2l1ELXmlUG+eQTuOceh4oeOnuIQ0mHuKqhaxTop0+Hl1+GZs2cX/e+x/dR886aVLuhmvMrN5QfXNHDLkBn77fffqNOnTosWrQIgKSkJKZNm3ZJuejoaL777js2bdqExWKhQ4cOdOzY8fxxi8XCunXr+OWXX3j55ZdZsmQJr7zyCtHR0Xz00UfO/04eIN+ehlLqGQAROaeUuvWiw6NdalVZ4+BBrTV+440OFT+WfIxtp7e5xJTMTJg8GR54wCXVU/226gQ3M6OXhhIi4vytANq0acOSJUsYN24cK1asIDQ076yYK1euZPDgwQQFBRESEsKgi+QTbr75ZgA6duxIbGysU34Kb6Og4anbc72/OEb0ehfYUnb57DMYPtyhEKX07HS61+vOA52df1ffuBE6dNCT4M5+kDv+v+Mcee8Iod1D8Q12fd5yg8GZNGvWjA0bNtCmTRuee+45XnnllTzLFaYKXqGCTiTm6+uLxZWqnx6kIKeh8nmf174hPywW+PxzLRviAJNXT+bdte863YyMDBg5Ep5/3vnzGAm/JBA7PpaqA1yf58NgcAXHjx8nODiYESNG8NRTT7Fx48Y8y/Xo0YOffvqJjIwMUlJSzg9nFURISAjJycnONtljFOQ0JJ/3ee0b8mPxYqhXD9q0caj4C1e9wP2d7ne6GePHQ2Qk3HGHc+tN3pDMrjt3EflDJMERZljKUDrZunUrXbp0oV27drz66qv897//zbNc586dueGGG4iKiuLmm2+mU6dO+Q5l5XD11VezY8eOMjMRnm8SJqWUFUhF9yqCgLScQ0CgiLhGCKkQSl0Spptu0nkzHJBBf3/t+1zZ8Eo61O7gdDP27NGpO6o5cX46PTadTVdsIuKjCKrfVN15FRvKHTt37qRly5aeNsMhUlJSqFSpEmlpaVx11VXMmDGDDh2c/z9bEvL6PZ2VhCnfgQoRMQPTJeXkSfjzTx3b6gAtq7ekVqVaTjXh3DkdXvvWW+DrxL+oJcXC1n5baTCugXEYhnLF2LFj2bFjBxkZGdx5551e5zBcjZE4dyUzZ8Itt0BI4aorG09s5JrG1+Dn49w/yRNPaAV2ZzoMAN+KvjR9vynh1zmeedBgKAvMmTPH0yZ4FMcTOhiKhohem+HAsJSIMPHviZxIPuFUE376CZYtg7ffdmq17HloDykbU4zDMBjKIaan4Sr+/FOH2Hbt6lDxH4b+4HQTDh/WnR0HOjpFovaY2gS3MpPeBkN5xPQ0XMWMGTB2bKELIlKyUug4oyMZFudlLBPRawkffBCuvNJp1XJsyjEOv3WYkA4h+AaaKS+DoTxinIYriI+HX3+FESMKLVopoBLf3fYdgX7O0yafNQvuvNO5YoTxP8Zz6NVDVB9iJr0NhvKMcRqu4Msv4YYboEqVAoudyzzHR+s+onEV56VlPHQInnxSm+Csye9z/5xj9927ab2wNUGNg5xTqcFQRsktblhSvvjiC44fP+6UupyFS52GUup6pdRupdQ+pdSzeRwfrpSKsW+rlVJRrrTHLYjA//6nh6YKITkzuVBZgqLy2mvw1FPQrp1z6kvfn862G7fR/PPmVO5U2TmVGgyG81gLGBIoV05DKeULTAH6Aa2AYUqpVhcVOwj0FJG2wERghqvscRsrV+p5jCuuKLBYWnYaIRVCeLjrw05r2maDDz7QPQ1nYEmxENMvhoYTGlJtoFGtNZRdYmNjadGiBXfffTetW7dm+PDhLFmyhCuuuIKIiAjWrVvHmTNnuPHGG2nbti3dunUjJiYGgISEBK677jrat2/Pvffee8GD4KxZs86vNL/33nvPO4hKlSoxfvx4unbtypo1a3jllVfo3LkzrVu3ZuzYsYgI8+fPJzo6muHDh9OuXTvS09PZsGEDPXv2pGPHjvTt25cTJ5wbcekQIuKSDegOLM61/xzwXAHlqwDHCqu3Y8eO4tWMGCHyzjuFFvtp908yZsEYpzW7ebNI9+4iVqvTqhSbzSaJfyc6r0KDIR927Nhxwf6E5RNkwvIJIiIS8UGE7I7fLdHHoqXD9A4iIvLEb0/I5FWTRUSk9uTacuzcMVl+cLn0/LyniIjcs/AemR49XUREKr1WSc5lnCuw/YMHD4qvr6/ExMSI1WqVDh06yOjRo8Vms8mCBQtk8ODB8tBDD8lLL70kIiJLly6VqKgoERF5+OGH5eWXXxYRkZ9//lkAiYuLkx07dsjAgQMlKytLRETuv/9+mTlzpoiIADJ37tzz7SckJJx/P2LECFm4cKGIiPTs2VPWr18vIiJZWVnSvXt3OX36tIiIfPPNNzJ69GiHfk97m9HijHu7MyrJs2IYAnySa38k8FEB5Z/KXf6iY2OBaCC6QYMGef5IXkFCgkhoqEhcnEPFrTbn3OHT00UiI0Xs16NT2DFyh5yLLvgfzWBwFnnd5NzJwYMHpWnTpuf3R44cKbNmzRIRkf3790tUVJS0a9dO9u/ff75MvXr15OzZsxIVFXXB51WqVJG4uDj58MMPpXbt2hIVFSVRUVHSrFkzmTBhgoiI+Pr6isViOX/O/PnzpUuXLtK6dWupU6eOvP766yJyodPYunWrhISEnK+vdevW0qdPnzy/jyudhivXaeQVa5rnAL5S6mrgLqBHXsdFZAb2oatOnTp5r1jirFnQv3+hAk/T1k8jPCicoa2HOqXZF16Ali21im1JEZuAggbPNSAowkx6G8oPObLmAD4+Puf3fXx8sFgs+OUhD52TjjmvtMwiwp133snrr79+ybHAwEB87ZEqGRkZPPDAA0RHR1O/fn1eeuklMjIuDcEXESIjI1mzZk3xvqCTcOVE+FGgfq79esAlMzpKqbbAJ8BgEUlwoT2uJWcC3IHsfP0i+tGlbhenNX3LLfDxx87JkbH/mf0cfe8oFVtWxMfPBNcZDDlcddVVzJ49G4A///yTatWqUbly5Qs+//XXX0lMTASgd+/ezJ8/n9OnTwNw5swZDh06dEm9OQ6iWrVqpKSkMH/+/PPHcsuqN2/enLi4uPNOIzs7m+3bt7vo2+aPK3sa64EIpVRj4Bg6qdMFwtxKqQbA98BIEdnjQltcz9q1OmlFr14FFltyYAlta7alRsUaJW7yzBmYOlX3NJzhMA6/dZgzv52h/Yr2Ja/MYChjvPTSS4wePZq2bdsSHBzMzJkzAZgwYQLDhg2jQ4cO9OzZkwYNGgDQqlUrJk2axHXXXYfNZsPf358pU6bQsGHDC+oNCwvjnnvuoU2bNjRq1IjOnTufPzZq1Cjuu+8+goKCWLNmDfPnz+eRRx4hKSkJi8XCY489RmRkpPt+BAqQRndK5Ur1B94DfIHPRORVpdR9ACLysVLqE+AWIMf9WqQQ6V6vlUYfM0aPET39dIHFJiyfwB1t7qB5teYlak4Ehg2DGjV0xFRJOTnzJAcnHKTDqg5UqFuh8BMMBidSmqTRSwMekUZ3BiLyC/DLRZ99nOv93UDhin7ezpkz8P33OmlFAaRmpfLy1S87pcnZsyEmBjZsKHldCYsS2D9uP+3+bGcchsFgKBAzaO0MPv8cBg3Sj/35kG3NpsOMDiSmJzqlye3b4euvIaiEc9UZRzLYNXoXrRe0pmKLik6xzWAwlF2Mym1Jsdn0xIJ9Iiw//H392XzvZoL8S3aXz8zUHZo8AjKKjDXNSmD9QDpu6EhgfedpXxkMhrKL6WmUlN9+g7CwAiXQT6ee5oFFDzhFlPDZZ2HSpBJXg81iY2P3jWSdzjIOw2AwOIzpaZSUKVO0BnkB4UtBfkEMbDYwz1juorBokZ462bSpRNWQdSoL/+r+tF/ZHr8QcwkYDAbHMT2NknDgAKxbp8OY8uFE8gliz8bSP6J/iZoSgZdf1qNg4SVImJdxJIMNXTeQtCLJOAyDwVBkjNMoCdOmwahRBc5Gb4/bzqK9i0rUjNUKFovWQuyR55p5x8g6lcWWa7dQ79F6hPUMK5FNBkN5JTY21iV5wj/++GO+/PJLp9frbMyjZnFJT9dRU+vW5Vsk05LJtU2u5dom15aoqZdfhuRkePfd4teRfSabLX22UHN4Teo/Xr/wEwwGQ57kOI077rij8MJF4L777nNqfa7C9DSKy5dfQvfu0KRJvkUeWPQAC3YtKFEzixZp3/TsJdlIHCc7MZuYvjGE9w2n4YsNCz/BYCiHXCxj/s8//9C2bVsyMjJITU0lMjKSbdu28eyzz7JixQratWvHu+++i9Vq5emnn6Zz5860bduW6dOnA1pqpFevXgwZMoQWLVowfPjw87Lpzz77LK1ataJt27Y89dRTgF5xPnnyZHbu3EmXLv/KDMXGxtK2bVuAsi2N7qrNK6TRrVaRZs1E/vyzwGLp2emSacksdjMJCSI1aoisWlXsKkRE5MSXJ2Tv43vFZrOVrCKDwUVcIo0+QWtw52zR0XrL/ZldMFZq1/73sw5aOV3uuefCsseOFd5+XjLmL7zwgjz55JPywAMPyGuvvSYiIsuXL5cBAwacP3f69OkyceJEERHJyMiQjh07yoEDB2T58uVSuXJlOXLkiFitVunWrZusWLFCEhISpFmzZuf/HxMTE+3feYK89dZbIiIXKOe+8cYbMnHiRK+RRjfDU8Xh55+hcmW46qo8D4sID/7yIK9c/QrVgoufvCg8HFatgqZNi3d+9plskjckU2tkLWSElDh6y2BwFy+9pLeLyUv1KK/EdjNm6M1Rli5dyoYNG87rPqWnp1OjRg3Gjx9P586dCQwM5IN89Hp+//13YmJizgsNJiUlsXfvXgICAujSpQv16tUDoF27dsTGxtKtWzcCAwO5++67GTBgAAMHDrykzttuu41vv/2WZ599lrlz5zJ37lx2797Ntm3b6NOnD6Az/tWuXdvxL+kkjNMoDpMn65yq+dyEBeGaxtcQHlS8MCcRuO8+uPVWuLYE0yGZRzM5t/Yc4X3CjcMwGApA8pExP3nyJCkpKWRnZ5ORkUHFipeqJogIH374IX379r3g8z///PMCuXVfX9/zEuvr1q1j6dKlfPPNN3z00UcsW7bsgnOHDh3Krbfeys0334xSioiICLZu3VrmpdHLJuvWweHDWo88D9Ky0/hj/x8MaTUEH1W8n3f6dFi9Wk+ZFIfshGwOvXqIiq0r0ujFRsWrxGAoR+QnYz527FgmTpzI8OHDGTduHHChXDlA3759mTZtGtnZ2QDs2bOH1NTUfNtKSUkhKSmJ/v37895777F58+ZLylx22WX4+voyceJEhg7VeXfKgzR62eTtt+GxxyCPhCwAh5MOszx2OX2b9s3zeGGsWwfjx+thqTweagol80SmnvTuF553GiyDwXAJecmYDx48GD8/P+644w6sViuXX345y5Yt48orr8TPz4+oqChGjRrFo48+SmxsLB06dEBEqF69OgsWLMi3reTkZAYPHkxGRgYiwrv5hEUOHTqUp59+moMHDwIQEBBQ9qXRXYFHpdEPHoTOnfVrSMglh89mnCUkIARfH99iN5GQADt2wJVXFv3c9APpbOmzhdpjatPg+QZmSMpQajDS6M7FldLoZniqKEyeDHffnafDAPi/Vf/Hp5s+LVbVKSkwejQEBBTPYaRsS2HTVZuo/2R9Gr7Q0DgMg8HgEszwlKMcOaK1yHftyrfIq9e8ilWsRa7aZoPhw3Vq8UqVim5aVpxe6d30nabUvKNm0SswGAwGBzE9DUd57TUYOzbPnBkWm4XrZ11PXFocfj5F98PPPQdJSVqVpKgdhKzTWQRUD6DD6g7GYRgMBpdjnIYjHDoE336rw2zzwM/Hj9d7v17svN9XXAHffaeHpoqCNc3Klj5bsJyzENSkhNmYDAaDwQHM8JQjvPqqXjhR7dKFeptObGLzyc2Mbj+6yNXOnQupqTq9eFEQm3Bqzilq3lGTjhs64uNnfL/BYHAPxmkUxsGDBeb/DvYPplalWkWudskSePhh/VoUrGlWdo7cSXZcNtUGVzPy5gaDwa2YR9TCmDQJHnggzyQWSw4soV7levSL6FekKjdv1ik45s0Duw6ZQ2SeyGRzz834VvQl6o8o4zAMhlLIqFGjzkuO5McXX3zB8bz0UbwA4zQKYutW+OknePzxPA//sPMH4tPii1SliBbGnT8fevZ0/Lyk1Uls7LKRqjdUpcXMFvhUMH86g6GsYpxGaUQEnnhCL8+uUuWCQ6lZqRxJOsKUAVNoGOa41PjGjdC7t17pXRSHkXE4g203byNiWgSNXmxk1mAYDC7gyy+/pG3btkRFRTFy5EgOHTpE7969adu2Lb179+bw4cOA7incf//9XH311TRp0oS//vqLMWPG0LJlS0aNGnW+vkqVKvHkk0/SoUMHevfuTVxc3CVt5iV1Pn/+fKKjoxk+fDjt2rUjPT3dOyTRc3CGVK47N7dJo//0k0iLFiJ2qeTcLNqzSB799dEiVbdhg0jNmiILFjh+Tva5bIn/OV5ERLISL7XDYCgrXCzlfWDCAVnO8ku2VbVXnT9+YMIBERFZVXtVnmVzH884llFg+9u2bZNmzZpJXFyciIgkJCTIwIED5YsvvhARkU8//VQGDx4sIiJ33nmnDB06VGw2myxYsEBCQkIkJiZGrFardOjQQTZt2iQiIoDMmjVLRERefvllefDBB8+fP2/evAKlznv27Cnr168XESmSJHp+v6fdHiON7jLS0uCRR+Djj8Hf/4JDZ9LP0D+iP/2aOj6PcfYsDByoqxs82LFzRARLkoXEpYmE9w/HP8y/8JMMhjJC45ca0/ilxgUez+Hy45cXWFdhxwGWLVvGkCFDqGaPkAwPD2fNmjV8//33AIwcOZJnnnnmfPlBgwahlKJNmzbUrFmTNm3aABAZGUlsbCzt2rXDx8fnvNjgiBEjuPnmmy9o01Gpc2+RRM/BOI28eOUV6NYNrrvugo+zrFn0+qIXf47602HZ8927oXlziInJM2L3EsQqHH7zMGm70mj5ZUuavlPMZBoGg8FhRArPN5P7eI7kuY+PzwXy5z4+PlgslkLPz2nTEalzR8u5CzOncTGbN8Nnn8E771zwcWpWKv4+/kSPjXbIYYjAhx/qOYyzZx1zGKk7UtnUYxOJfyTS+NX8n7IMBoNz6d27N99++y0JCQmAlka//PLL+eabbwCYPXs2PXr0KFKdNpvtfJTUnDlzLjm/IKnz3PLr3iKJnoPpaeQmPV2LQL37LtS6cO3F80ufp2u9rtzRpvBk8hkZOkp33TpYsQLCwgoub8uycfiNwxz78BiNJjaiztg6KB8z2W0wuIvIyEheeOEFevbsia+vL+3bt+eDDz5gzJgxvPXWW1SvXp3PP/+8SHVWrFiR7du307FjR0JDQ5k7d+4FxwuSOh81ahT33XcfQUFBrFmzxisk0XMw0ui5eeQRiIuDOXPOi0BZbVbSstPw8/EjwDegUNlzEcjM1IvIx41zTIBw87Wb8Qn0odm0ZgTWD3TGNzEYShVlURq9UqVKpKSkeKRtI43uDr7+Wuf+njr1AtXAeTvm8cKyFwjyDyrUYfz4I/ToAb6+MHFiwQ4j+2w2+8ftR0Ro8UUL2vzUxjgMg8Hg9ZjhKdDzGI88ojU97GsyLDYLsWdjGRo5lBtb3Fjg6WfOwKOPwpo18PnnlwRcXUDGoQwyjmQQenkowc2DEasQWM84C4OhrOGpXoarMT2NvXthwACtSx4Vdf7j6OPRvPLXKyilCPTL+6YuAllZsH+/9jVbtuSfQCn9YDq779tNdIdoUjakoHwUtcfUNmKDBoOd0jZU7q24+ncs3z2N/fuhTx8dYjtkCAD7zuxjw/ENDG09lK51u+Z76qpV8PTTcMst8OSTOgvsxYgIZ5ed5egHR0lalUSdsXXosrsLAdWKqIFuMJRxAgMDSUhIoGrVqkbxoASICAkJCQQGum70ovw6jbVr4aab4KWX4K67sImNTEsmCkVKlu5W5nXxisBtt+nTJ02CESMurTr9QDr+Vf1J2ZrCvsf2UffhurSa0wrfisXPHW4wlGXq1avH0aNH85TaMBSNwMBA6tWr57L6y1/0lNWqc31PngwzZ0L//gC8sfINfJUvT1/x9CWnJCVpgcGTJ+GFF2DxYrjqKgiy5z0Sm5CyJYX0fenUuLUGex/eS/XbqhPaIxTI2/kYDAaDOykV0VNKqeuVUruVUvuUUs/mcVwppT6wH49RSnVwmTEi8Msvehzp119h/Xro359b593K4aTDPNr1UZ68/MnzRXPmsEaOhPr1YeFC6GofrerTy0bWliSOvn+UHcN2sLrOanYM3UHKZn1SxIcRhF0ZhlLKOAyDwVCmcNnwlFLKF5gC9AGOAuuVUgtFZEeuYv2ACPvWFZhmf3UOycmwYQMsWwZff80/dWx0fuFNlrUL5bttb/JhvWncFfECFbLqEFjZj0mTtC/Z9I+VAd2zePOZLB5oncXbjwRTvWMw22/bjq1XK07NOcWxD45RuVtlqvStQuNJjQm6zKRbNRgMZR9Xzml0AfaJyAEApdQ3wGAgt9MYDHxpV2Bcq5QKU0rVFpECdX9jnnqfZVPbEmLxI8vHRkiFBNr4JrI6sy7LfENJErjGlsSz3MJ/a7xJg+PXEeB3A2lHssic5IdvtqLergcIBO4IrMa5W05z28w6dHp3LT3OZaJ8ocKmAPY9EUClWgH4tqmD8qlI7dG1wQa1R9fW7w0Gg6Gc4bI5DaXUEOB6Ebnbvj8S6CoiD+Uq8zPwhoistO8vBcaJSPRFdY0Fxtp3WwPbXGK0c6kGFC1Dk2cwdjqX0mBnabARjJ3OprmIhJS0Elf2NPIazL/YQzlSBhGZAcwAUEpFO2Myx9UYO52LsdN5lAYbwdjpbJRSTtFfcuVE+FGgfq79esDF+QsdKWMwGAwGL8GVTmM9EKGUaqyUCgBuBxZeVGYh8B97FFU3IKmw+QyDwWAweA6XDU+JiEUp9RCwGPAFPhOR7Uqp++zHPwZ+AfoD+4A0YLQDVc9wkcnOxtjpXIydzqM02AjGTmfjFDtL3eI+g8FgMHgOo5ZnMBgMBocxTsNgMBgMDuNVTqMksiOFnetGG4fbbYtRSq1WSkXlOharlNqqlNrsrPC3EtjZSymVZLdls1JqvKPnutnOp3PZuE0pZVVKhduPufP3/EwpdVoplecaIS+5Nguz0VuuzcLs9JZrszA7veXarK+UWq6U2qmU2q6UejSPMs67PkXEKzb0ZPl+oAkQAGwBWl1Upj/wK3p9RzfgH0fPdaONlwNV7O/75dho348FqnnJb9kL+Lk457rTzovKDwKWufv3tLd1FdAB2JbPcY9emw7a6PFr00E7PX5tOmKnF12btYEO9vchwB5X3ju9qadxXnZERLKAHNmR3JyXHRGRtUCYUqq2g+e6xUYRWS0iifbdtei1J+6mJL+Hu37L4rQ1DPjaRbYUiIj8DZwpoIinr81CbfSSa9OR3zI/3HltFtVOT16bJ0Rko/19MrATqHtRMaddn97kNOoCR3LtH+XSL55fGUfOdZeNubkL7d1zEOB3pdQGpaVRXIWjdnZXSm1RSv2qlIos4rnOwOG2lFLBwPXAd7k+dtfv6QievjaLiqeuTUfx9LXpMN50bSqlGgHtgX8uOuS069ObkjCVRHbEITkSJ+BwO0qpq9H/mD1yfXyFiBxXStUA/lBK7bI/zXjCzo1AQxFJUUr1Bxag1Ybd9VtSxLYGAatEJPeTn7t+T0fw9LXpMB6+Nh3BG67NouAV16ZSqhLacT0mIucuPpzHKcW6Pr2pp1ES2RF3yZE41I5Sqi3wCTBYRBJyPheR4/bX08AP6K6hKyjUThE5JyIp9ve/AP5KqWqOnOtOO3NxOxd1/934ezqCp69Nh/CCa7NQvOTaLAoevzaVUv5ohzFbRL7Po4jzrk93TNQ4OJnjBxwAGvPvhEzkRWUGcOFkzjpHz3WjjQ3QK9wvv+jzikBIrver0SrAnvota/Hv4s4uwGH77+qW37IofzcgFD22XNETv2euNhuR/+StR69NB230+LXpoJ0evzYdsdNbrk37b/Ml8F4BZZx2fXrN8JSUQHYkv3M9ZON4oCowVemsfRbRCpg1gR/sn/kBc0TkN2fbWAQ7hwD3K6UsQDpwu+iryC2/ZRHsBLgJ+F1EUnOd7rbfE0Ap9TU6qqeaUuooMAHwz2WnR69NB230+LXpoJ0evzYdtBO84NoErgBGAluVUpvtnz2Pfkhw+vVpZEQMBoPB4DDeNKdhMBgMBi/HOA2DwWAwOIxxGgaDwWBwGOM0DAaDweAwxmkYDAaDwWGM0zAY7NhVSnMUS+fZ5SFQStVUSs1RSh2wy0KsUUrddNG57yuljimlfC76/MbcKq35tFtdKeXKkEyDwWkYp2Ew/Eu6iLQTkdZAFnCf0sH2C4C/RaSJiHRErwA+L/ZndxQ3oTV8rrqozmeAqQU1KiJxwAml1BVO+yYGg4swTsNgyJsVQFPgGiAr12IuROSQiHyYq+zVwDZgGlrtFAClVDMgU0Ti7fuXKaXWKqXWK6VeUUql5KpjATDcZd/GYHASxmkYDBehlPJD55vYCkSiBfQKIkcW+wdgoF0HCPRK3dznvg+8LyKduVTfJxq4soSmGwwuxzgNg+FfguwyDNFovaNPLy6glJpil+xeb98PQMszLBCtLPoPcJ29eG0gLtfp3YF59vdzLqr6NFDHSd/DYHAZXqM9ZTB4Aeki0i73B0qp7cAtOfsi8qBdcTUnhef1aNG6rXatoWC0ts8itG5SqINtB9rLGwxejelpGAwFswwIVErdn+uz4FzvhwF3i0gjEWmEVgu9zh55tRM9L5LDWv51QLdf1E4z9LyIweDVGKdhMBSAXV31RqCnUuqgUmodMBMYZ3cMfdG9ipzyqcBKdGKev4H29ggsgMeAJ+x11AaScjV1de56DAZvxajcGgwuRCn1PvCTiCyxO5l0ERGl1O3AMBEZbC/3NzoxUmJB9RkMnsbMaRgMruU1oKv9fUfgI3vP4ywwBvTiPuAd4zAMpQHT0zAYDAaDw5g5DYPBYDA4jHEaBoPBYHAY4zQMBoPB4DDGaRgMBoPBYYzTMBgMBoPD/D8EPX5MtFoL4QAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 绘制易损性曲线\n",
    "# 确定坐标轴\n",
    "plt.xlim((0,2))\n",
    "plt.ylim((0,1))\n",
    "\n",
    "plt.title(\"Fragility of Bridge\") \n",
    "plt.xlabel(\"PGA(g)\") \n",
    "plt.ylabel(\"Exceeding Probability\")\n",
    "\n",
    "# slight\n",
    "x1 = np.linspace(0.01,2,200)\n",
    "y1 = []\n",
    "for i in x1:\n",
    "    k = stats.norm.cdf((math.log(i) - math.log(slightMedian))/(0.6))\n",
    "    y1.append(k)\n",
    "\n",
    "plt.plot(x1,y1,'r-',linewidth=1,linestyle='-',label='slight')\n",
    "\n",
    "# moderate\n",
    "x2 = np.linspace(0.01,2,200)\n",
    "y2 = []\n",
    "for i in x2:\n",
    "    k = stats.norm.cdf((math.log(i) - math.log(moderateMedian))/(0.6))\n",
    "    y2.append(k)\n",
    "\n",
    "plt.plot(x2,y2,'g-',linewidth=1,linestyle=':',label='moderate')\n",
    "\n",
    "# extensive\n",
    "x3 = np.linspace(0.01,2,200)\n",
    "y3 = []\n",
    "for i in x3:\n",
    "    k = stats.norm.cdf((math.log(i) - math.log(extensiveMedian))/(0.6))\n",
    "    y3.append(k)\n",
    "\n",
    "plt.plot(x3,y3,'b-',linewidth=1,linestyle='--',label='extensive')\n",
    "\n",
    "# complete\n",
    "x4 = np.linspace(0.01,2,200)\n",
    "y4 = []\n",
    "for i in x4:\n",
    "    k = stats.norm.cdf((math.log(i) - math.log(completeMedian))/(0.6))\n",
    "    y4.append(k)\n",
    "\n",
    "plt.plot(x4,y4,'m-',linewidth=1,linestyle='-.',label='complete')\n",
    "\n",
    "# 图例和绘制\n",
    "plt.legend(['slight','moderate','extensive','complete'],loc = 'lower right')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 计算超越概率\n",
    "exceedingProbability = []\n",
    "exceedingProbability.append(stats.norm.cdf((math.log(PGA) - math.log(slightMedian))/(0.6)))\n",
    "exceedingProbability.append(stats.norm.cdf((math.log(PGA) - math.log(moderateMedian))/(0.6)))\n",
    "exceedingProbability.append(stats.norm.cdf((math.log(PGA) - math.log(extensiveMedian))/(0.6)))\n",
    "exceedingProbability.append(stats.norm.cdf((math.log(PGA) - math.log(completeMedian))/(0.6)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 计算处于五种破坏状态的概率\n",
    "damageState = [0]*5\n",
    "damageState[0] = exceedingProbability[3]\n",
    "damageState[1] = exceedingProbability[2] - exceedingProbability[3]\n",
    "damageState[2] = exceedingProbability[1] - exceedingProbability[2]\n",
    "damageState[3] = exceedingProbability[0] - exceedingProbability[1]\n",
    "damageState[4] = 1 - exceedingProbability[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "无损伤的概率是: 0.3549814490622353\n",
      "轻微损伤的概率是: 0.14501855093776472\n",
      "中等损伤的概率是: 0.21252864582899444\n",
      "严重损伤的概率是: 0.12384052434297221\n",
      "完全损伤的概率是: 0.16363082982803334\n"
     ]
    }
   ],
   "source": [
    "# 从完全毁坏到无破坏输出概率列表\n",
    "print('无损伤的概率是:',damageState[4])\n",
    "print('轻微损伤的概率是:',damageState[3])\n",
    "print('中等损伤的概率是:',damageState[2])\n",
    "print('严重损伤的概率是:',damageState[1])\n",
    "print('完全损伤的概率是:',damageState[0])"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3.9.9 ('pyincoreEnv')",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.9"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "65f8545d5b179584e4344a378342d3091d162296694f191a9d66db391bbacd9c"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
